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ABSTRACT 


The  non-linear  theory  of  one-dimensional  consolidation 
formulated  by  Davis  and  Raymond  (1965)  corresponds  with  experi¬ 
mental  results  for  normally  consolidated  soils. 

This  thesis  presents  an  analytical  study  of  two-dimensional 
consolidation  of  normally  consolidated  soils.  The  assumptions 
adopted  in  the  Davis  and  Raymond  theory  (1965)  and  the  work  of 
Henkel  and  Sowa  (1963)  on  normally  consolidated  soils  are  used  in 
this  thesis  to  formulate  general  equations.  The  total  stress  com¬ 
ponents  are  derived  from  elastic  stress  distribution  for  a  strip 
load  on  a  half-space. 

A  general  numerical  method  of  solution  is  presented  for 
thin  layers  in  which  total  stress  remains  constant  with  respect 
to  time  and  space. 

Discussions  are  presented  on  the  significance  of  non¬ 
elastic  effects  on  the  rate  of  dissipation  of  excess  pore  pressure. 
These  discussions  indicate  that  the  rate  of  consolidation  considered 
is  almost  independent  of  the  pore  pressure  induced  by  shear  stresses. 
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CHAPTER  I 

INTRODUCTION 

The  theory  of  consolidation,  as  formulated  by  Terzaghi 
(1925),  is  based  on  the  development  and  subsequent  dissipation  of 
excess  pore  pressure  in  a  saturated  soil  medium.  Consolidation  is 
a  process  involving  a  time  dependent  volume  reduction  and  associated 
decrease  in  moisture  content.  In  deriving  the  differential  equation 
governing  one-dimensional  consolidation,  Terzaghi  assumed  constant 
coefficients  of  permeability  and  compressibility  for  the  soil  and 
small  strains.  All  these  simplifying  assumptions  make  possible 
a  linear  theory  and  the  governing  equation  facilitates  closed  form 
solution  for  given  boundary  conditions. 

One-dimensional  consolidation  is  due  either  to  a  lateral 
physical  restraint,  as  in  the  oedometer  test,  or  due  to  a  surface 
loading  which  is  of  uniform  magnitude  and  extends  infinitely  in  all 
directions.  Published  field  evidence  of  the  rate  at  which  foundations 
on  clay  settle  suggests,  that  the  actual  rates  are  generally  faster 
than  those  predicted  by  the  one-dimensional  consolidation  theory. 

It  is  obvious  that  in  many  practical  cases,  the  geometric  conditions 
are  not  one-dimensional  and  that  horizontal  dissipation  of  pore 
pressure,  which  is  not  taken  into  account  in  the  one-dimensional 
theory,  increases  rate  of  settlement.  In  spite  of  these  limitations, 
the  mathemati cal  simplicity  of  one-dimensional  theory  permits  a 
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large  variety  of  solutions  to  be  obtained.  The  use  of  three- 
dimensional  theory  for  practical  foundation  problems  appears  to  be 
hindered  by  a  shortage  of  theoretical  solutions  covering  an  adequate 
range  of  practical  situations. 

The  world  is  non-linear  and  contrary  to  many  notions  in¬ 
herited  from  classical  mechanics,  the  events  that  take  place  seldom 
are  linear.  There  are  no  strictly  homogeneous  soils.  All  defor¬ 
mations,  velocities  and  stress  increments  need  not  be  small.  There 
are  no  ideal  fluids.  This  does  not  imply  that  all  of  these  com¬ 
plexities  of  nature  need  be  considered  in  practical  designs,  but 
the  best  design  will  always  emerge  from  the  most  accurate  assess¬ 
ments  of  natural  condition.  Terzaghi's  linear  theory  relies  on  a 
set  of  unreal  assumptions  and  is  therefore  only  an  approximation. 

Non-linear  behaviour  in  soil  systems  is  usually  placed 
into  one  of  three  categories. 

1)  Geometric  non-linearity,  which  arises  from  non-linear 
terms  in  the  kinematic  equations, 

2)  Material  non-linearity ,  which  arises  from  non-linearities 
in  the  constitutive  equations, 

3)  Combined  geometric  and  material  non-linearity. 

Differences  in  observed  and  predicted  results  have  stimu¬ 
lated  many  workers  to  introduce  a  series  of  modifications  to  Terzaghi's 
theory.  Schiffman  and  Gibson  (1964)  studied  the  progress  of  con¬ 
solidation  in  a  stratum  of  clay  through  which  there  is  a  known 
continuous  variation  of  the  consolidation  parameters.  They  have 
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presented  a  theory  taking  these  changes  into  account.  To  overcome 
the  inaccuracies  incorporated  in  the  Terzaghi  theory,  Davis  and 
Raymond  (1965)  presented  a  small  strain  theory  for  compressible 
normally  consolidated  clay  by  taking  into  consideration  a  stress 
dependent  coefficient  of  compressibility  and  coefficient  of  perme¬ 
ability.  Barden  and  Berry  (1965)  also  presented  a  small  strain 
theory  with  a  variation  of  consolidation  parameters.  Large  strains 
have  been  considered  by  Mikasa  (1965).  Gibson  et  al  (1967)  developed 
a  general  theory  of  one-dimensional  consolidation  for  which  the  limit¬ 
ation  of  small  strains  is  not  imposed  and  both  the  variation  of  soil 
compressibility  and  permeability  during  consolidation  has  been  taken 
into  account. 

Unavoidably,  many  of  the  above  modifications  lead  to  non¬ 
linear  partial  differential  equations  which  render  classical  methods 
of  mathematical  analysis  inapplicable.  When  such  is  the  case,  re¬ 
course  is  usually  made  to  numerical  methods  for  quantitative  solu¬ 
tions  to  problems. 

1.1  Terzaghi 's  Classical  One-Dimensional  Theory 

For  many  field  problems,  the  one-dimensional  consolidation 
theory  (Terzaghi,  1925)  is  sufficiently  accurate.  The  basic  as¬ 
sumptions  involved  in  the  theory  are: 

(i)  the  clay  is  saturated  and  homogeneous 

(ii)  the  pore  water  and  soil  grains  are  incompressible  in 
comparison  with  the  soil  skeleton 

(iii)  the  flow  and  consolidation  are  one-dimensional 
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(iv)  the  fluid  flow  follows  Darcy's  Law 

(v)  the  coefficients  of  permeability  and  compressibility 
are  constant  during  the  process  of  consolidation 

(vi)  strains  in  the  soil  skeleton  are  controlled  exclusively 
by  the  effective  stress  via  a  linear  time  independent 
relation 

(vii)  strains,  velocities  and  stress  increments  are  small. 

The  basic  theory  considers  only  the  diffusion  of  water 
through  a  porous  medium.  Equating  the  volume  of  pore  water  squeezed 
out  with  the  reduction  in  void  ratio  of  an  element  of  the  soil,  the 
well  known  linear  equation  of  Terzaghi  one-dimensional  consolidation 
is  deduced.  The  governing  equation  may  be  written  as 

2 

3u  _  r  3  u 
3t  v3z2 

where  Cv  denotes  coefficient  of  consolidation  and 
u  denotes  pore  water  pressure. 

The  application  of  one-dimensional  theory,  which  relies 
on  a  set  of  unreal  assumptions  concerning  material  properties  and 
the  dimensionality  of  compression,  paradoxically,  yields  acceptable 
results  and  permits  a  wide  range  of  solutions  to  be  obtained.  An 
examination  of  Equation  1.1  reveals  that  it  is  linear  and  closed 

i 

form  solutions  are  possible.  An  important  feature  of  the  linearity 
is,  the  additive,  or  superposable ,  nature  of  solutions.  Although 
Equation  1.1  describes  the  behaviour  of  the  soil  in  an  artificial 
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way,  its  superposition  properties  and  closed  form  solutions  render 
it  both  attractive  and  useful,  provided  the  inherent  limitations 
are  observed. 

Because  of  linearity  it  is  appropriate  to  normalize  the 
various  parameters  to  make  the  equation  non-dimensional.  In  this 
way,  the  solution  obtained  in  terms  of  dimensionless  parameters  is 
in  a  more  suitable  form  for  general  application.  Recent  develop¬ 
ments  in  numerical  analysis  and  the  modern  advances  in  computer 
technology  provide  improved  capability  of  achieving  greater  versatility 
in  handling  one-dimensional  consolidation  equations. 

1.2  Terzaghi-Rendul ic  Two-Dimensional  Theory 

When  a  soil  consolidates  under  an  applied  load,  compression 
and  drainage  often  take  place  in  three  dimensions  rather  than  in  one 
as  assumed  in  Terzaghi's  one-dimensional  theory.  However,  the  exist¬ 
ing  knowledge  of  the  consolidation  characteristics  of  soil  in  three- 
dimensions  is  still  in  the  preliminary  stages.  Since  the  displace¬ 
ment  of  principal  importance  in  structural  applications  is  vertical, 
with  due  reservations,  the  influence  of  other  than  vertical  directions 
may  often  be  neglected.  However,  the  process  of  drainage  of  water 
in  three  dimensions  does  have  an  effect  on  the  rate  of  settlement 
of  a  soil  loaded  at  the  surface. 

,,  The  simple  consolidation  theory  generally  attributed  to 

Terzaghi  (1925)  and  Rendulic  (1937)  considers  the  three  dimensional 
diffusion  of  water  through  a  porous  medium.  The  conventional  theory 
of  consolidation,  as  formulated  by  Terzaghi  (1925),  assumes  that 


6 


only  the  excess  pore  pressure  contributes  to  the  progress  of  con¬ 
solidation.  For  a  time-independent  loading  the  progress  of  con¬ 
solidation  is  governed  by  a  single  relationship  of  heat  conduction 
type  as  Equation  1.1. 

In  conventional  theory  the  internal  total  stresses  are 
generated  only  by  the  applied  loads.  This  interpretation  has  led 
to  the  Terzaghi -Rendul ic  (pseudo  three-dimensional)  theory.  For  a 
time-independent  loading,  the  excess  pore  pressures  are  governed  by 

CPV2U=|^  1.2 

where  Cp  denotes  coefficient  of  consolidation. 

The  formulation  of  pseudo  three-dimensional  theory  involves 
theoretical  inconsistencies  as  it  takes  no  account  of  strain  com- 
patability  within  the  soil  mass.  As  with  one-dimensional  theory, 
during  the  process  of  consolidation,  the  total  stress  which  may  be 
calculated  from  elastic  stress  distribution  theory,  remains  constant 
for  time-independent  loading. 

Equation  1.2  is  a  half-way  point  between  the  one-dimensional 
theory  and  a  true  three-dimensional  theory.  This  equation  retains 
all  the  characteristics  of  one-dimensional  theory  and  is  based  on 
the  same  assumptions.  Equation  1.2  provides  a  means  for  dissipation 
of  excess  pore  pressure  in  all  directions  but  the  consolidation  is 
restricted  to  the  vertical  direction  only. 

Under  two-dimensional  drainage  conditions  Equation  1.2 
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becomes 


fd^u  ,  9^ux  9u 

+  —7>  3t 


9x 


9z 


1.3 


Equation  1.3  governing  the  rate  of  dissipation  is  of  the 
heat  conduction  type  and  may  be  used  for  all  types  of  problems 
solved  using  the  classical  one-dimensional  theory.  Since  Equation 
1.3  has  all  the  characteristics  of  the  classical  one-dimensional 
theory,  it  is  called  the  Terzaghi  equation  for  consolidation  in  two- 
dimensi ons . 

1.3  Non-Linear  Consolidation  Theories 

To  improve  correspondence  between  prediction  and  performance, 
many  of  the  assumptions  employed  in  Terzaghi's  theory  have  been 
modified  by  various  workers  from  time  to  time.  Modest  deviations 
from  the  basic  theory  yield  non-linear  partial  differential  equations 
which  render  classical  methods  of  mathematical  analysis  inapplicable. 
Complexities  of  non-linear  analysis  have  created  a  "non-linear  barrier" 
which,  only  in  very  recent  times,  has  begun  to  be  penetrated  by  means 
of  the  efficient  use  of  modern  high  speed  computers.  Faced  with 
the  necessity  of  obtaining  solutions,  when  classical  methods  show 
sign  of  breakdown,  recourse  is  usually  made  to  numerical  techniques. 

Finite  Strain  Non-Linear  Theory  due  to  Gibson  et  al  (1967) 

Gibson  et  al  (1967)  published  a  general  theory  of  one¬ 
dimensional  consolidation  for  which  the  limitations  of  small  strains 
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are  not  imposed,  and  the  variation  of  compressibility  and  perme¬ 
ability  during  consolidation  has  been  taken  into  account.  The 
validity  of  Darcy's  Law  has  been  assumed  but  it  has  been  recast 
in  a  form  in  which  the  relative  velocity  of  soil  skeleton  and  the 
pore  fluid  related  to  the  excess  pore  pressure  gradient. 

For  the  particular  case  of  a  horizontal  thin  layer,  the 
self  weight  effect  is  negligible  and  the  governing  equation  in 
terms  of  the  void  ratio,  e,  is 


9  r  p  9e  ~i  _  9e 
9a  llF  3aJ  "  9t 


1.4 


where  "a"  is  the  vertical  distance  of  an  element  of  the  soil  skeleton 
from  an  embedded  datum  plane  (Fig.  1.1).  The  quantity  Cp  is  closely 
related  to  the  familiar  coefficient  of  consolidation  Cy  of  Terzaghi's 
theory,  by  the  expression 
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By  assuming  Cp  to  remain  constant  during  consolidation, 
the  differential  Equation  1.4  may  be  linearized  under  a  specific 
pressure  increment,  so  that 


9e 

3t 


1.6 


The  important  feature  of  Equation  1.6  is  that,  unlike 
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the  corresponding  Terzaghi  equation,  it  holds  without  any  restriction 
on  the  form  of  relation  between  effective  pressure  and  void  ratio. 

The  linearization  procedure  may  lead  to  errors,  either 
because  of  the  nature  of  the  soil  or  the  size  of  the  load  increment. 
Equation  1.4  must  be  solved  numerically  using  a  relation  between 
Cp  and  e  determined  experimentally. 

Davis  and  Raymond's  Non-Linear  Theory 

It  is  known  that  the  Terzaghi  theory  does  not  hold  for 
the  case  of  soft  clay  whose  permeability,  k,  and  compressibility,  m  , 
change  during  consolidation.  Davis  and  Raymond  (1965)  introduced  a 
modification  to  the  Terzaghi  theory  in  order  to  overcome  the  in¬ 
accuracies  of  the  following  assumptions: 

(i)  the  coefficient  of  permeability,  k,  is  a  constant 

during  consolidation  under  a  given  stress  increment;  and 
(ii)  the  coefficient  of  compressibility,  m  ,  remains  constant 
during  consolidation  under  a  given  stress  increment. 

Instead  of  assuming  that  m  is  constant,  they  adopt  the 
well  established  empirical  linear  relationship  between  void  ratio, 
e  and  the  logarithm  of  the  effective  pressure,  o'  as 


where  e  denotes  the  void  ratio  of  soil  subjected  to  o' 
o  0 

a'  denotes  the  effective  pressure  at  point  "o"  on  the 
e  -  log  o'  curve 
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and  Ic  denotes  the  compression  index  of  a  soil  and  is  constant. 

The  coefficient  of  compressibility,  mv  of  the  soil  skeleton  is  given 
by 


-  1  8e  ,  0 

v  '  (1+e)  8a 1  1,8 

Differentiating  Equation  1.7  with  respect  to  effective  pressure, 
a',  and  substituting  in  Equation  1.8,  yields 

0.434  I 

mv  =  (1+e)  o'  1,9 

During  the  consolidation  process  the  total  volume  (1+e)  varies  less 
with  time  than  the  effective  pressure,  a',  hence  (1+e)  may  be  con¬ 
sidered  constant  for  any  load  increment.  With  this  assumption 
Equation  1.9  becomes 


1.10 


where  A  is  a  constant  which  depends  on  the  initial  void  ratio  and 
compressibility  index,  Ic. 

Davis  and  Raymond  point  out  that  during  any  given  pressure 
increment  the  coefficient  of  consolidation,  Cy,  usually  changes  much 
less  than  the  coefficient  of  compressibility,  m  ,  for  normally  con¬ 
solidated  clay.  They  assumed  that  Cy  remains  constant,  i.e., 
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CV  =  =  constant  1.11 

\r  w 

where  y  denotes  unit  weight  of  water, 
w  v 

This  can  be  seen  to  be  equivalent  to  assuming  that  as  the  soil  parti¬ 
cles  move  closer  together  the  decrease  in  permeability,  k,  is  pro¬ 
portional  to  the  decrease  in  compressibility,  m  . 

By  adopting  Equations  1.10  and  1.11  in  place  of  the  as¬ 
sumptions  (i)  and  (ii)  in  the  Terzaghi  theory,  Davis  and  Raymond 
derived  the  following  governing  equation  which  is 

2 

1  3a 1  _  p  rl3Up  1  3u  3a 1 -i  n  i0 

a'  3t  ”  "  LV  La'  -  2  u  ,  n2  3z  3z  J  l-u 

a  z  ^a  ) 

Equation  1.12  is  a  weakly  non-linear  equation  for  a  normally  con¬ 
solidated  soil  under  a  constant  load.  Using  the  substitution, 

i  c'^-u 

w  =  log10  a7!  =  log10  "aT7~  1,13 

t  t 

where  a'f  denotes  final  effective  pressure,  and  differentiating 
this  expression  with  respect  to  z  and  t  and  substituting  in  Equation 
1.12  yields  a  simpler  version  of  the  differential  equation  in  terms 
of  the  function,  w 

r  32w  _  3w 

v  az2  "  at 


1.14 
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This  equation  is  identical  in  form  to  that  of  the  Terzaghi 
linear  theory  and  a  solution  to  the  equation  will  be  unique  depend¬ 
ing  on  the  specified  boundary  and  initial  conditions. 

Modification  by  Barden  and  Berry  (1965) 

Barden  and  Berry  (1965)  modified  the  assumptions  in 
Terzaghi's  theory  to  agree  more  closely  with  the  physical  properties 
of  such  materials  as  clays  which  exhibit  marked  secondary  consoli¬ 
dation  (creep)  effects.  The  modifications  are: 

(i)  the  effects  of  structural  viscosity  which  controls  creep 
require  that  void  ratio  be  a  function  of  both  effective 
stress  and  time, 

(ii)  for  normally  consolidated  clays,  the  void  ratio-logarithm 
of  effective  stress  relation  is  a  better  approximation, 

(iii)  a  void  ratio-logarithm  of  permeability  relation  is  con¬ 
sidered. 

Combining  e-logk  relation  with  the  e-loga'  relation  a  general  re¬ 
lation  for  permeability  may  be  developed 

k  =  kf  (1  +  bun)  1.15 

where  kf  denotes  final  permeability 
and  b  and  n  are  constants. 

Satisfying  the  equation  of  continuity  for  one-dimensional 
vertical  flow  in  a  saturated  soil,  Barden  and  Berry  have  presented 


< 


. 
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a  small  strain  theory  with  the  following  non-linear  equation 


JL 

9z 


*f 


'Y, 


w 


(1  +  bu11 ) 


9  u  -j 
9YJ 


CC  _J_  9_u 

(l+e0)  (a'f-uy  9t 


1.16 


Exact  solutions  to  Equation  1.16  are  difficult  to  obtain, 
hence  finite  difference  solutions  have  been  developed. 


1.4  Scope  of  the  Study 

The  non-linear  theory  of  Davis  and  Raymond  (1965)  is  in 
excellent  agreement  with  experimental  results  for  large  as  well  as 
small  pressure  increment  ratios  (Burland  and  Roscoe,  1969).  In 
engineering  practice  plane  flow  of  pore  water  is  assumed  when 
analysing  dissipation  of  pore  water,  for  example,  in  rolled  filled 
layers  of  earth  dams  or  strip  foundations,  or  embankments.  An 
analytical  approach  is  made  to  study  the  two-dimensional  consoli¬ 
dation  of  normally  consolidated  soils  (a  strip  load  of  width  2B  of 
uniform  intensity  q  resting  on  the  half-plane).  Parameters  for 
normally  consolidated  soil  considered  herein  are  likely  to  arise 
in  practice. 
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CHAPTER  II 

TWO-DIMENSIONAL  CONSOLIDATION  OF  NORMALLY 
CONSOLIDATED  SOILS 
GENERAL  APPROACH 

2.1  Method  of  Analysis 

Some  of  the  well  known  and  generally  accepted  principles 
of  soil  mechanics  and  foundation  engineering  are  employed  in  this 
study.  Classical  theory  of  elasticity  is  used  which  assumes  that 
a  linear  relation  holds  between  stress  and  strain  and  that  strains 
are  small.  Numerical  techniques  are  utilized  in  conjunction  with 
the  assumptions  used  in  the  one-dimensional  Davis  and  Raymond  non¬ 
linear  theory. 

2.2  Assumptions  Employed 

The  following  physical,  mathematical  and  engineering  as¬ 
sumptions  are  made  in  the  study. 

(a)  Normally  consolidated,  completely  saturated  soil. 

(b)  Darcy's  law  is  valid. 

(c)  Problem  is  one  of  plane  strain. 

(d)  Bishop's  (1966)  relationship  is  used  to  evaluate  the 
intermediate  principal  stress,  a 

(e)  Load  increment  is  applied  instantaneously. 
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(f)  Initial  excess  pore  pressures  can  be  expressed  in  terms 
of  octahedral  normal  stresses,  octahedral  shear  stresses 
and  the  Skempton  pore  pressure  parameter,  A. 

(g)  Volume  change  is  mainly  due  to  hydrostatic  stress. 

(h)  Empirical  relationship  between  void  ratio,  e,  and 
logarithm  of  the  average  effective  stress,  0',  (Henkel 
and  Sowa,  1963)  is  valid. 

(i)  In  spite  of  the  variation  of  the  coefficient  of  com¬ 
pressibility  and  permeability  with  the  changes  in  stress, 
their  combined  effect  is  balanced  to  maintain  the  co¬ 
efficient  of  consolidation  constant. 

(j)  Soils  are  homogeneous. 

(k)  The  coefficient  of  lateral  earth  pressure  is  constant 
and  equal  to  the  "at  rest"  value  (k  ). 

(l)  Representative  values  for  stresses  are  obtained  for  a 
strip  load  using  the  Bousinesq  relations. 

(m)  The  coefficient  of  consolidation  is  independent  of  the 
geometry  of  the  loading  and  the  boundary  conditions. 

(n)  During  the  consolidation  process,  the  load  distribution 
on  the  surface  of  the  soil  layer  and  the  associated 
total  stress  distribution  in  the  soil  mass  are  unaffected. 

(o)  The  ground  water  table  is  at  the  upper  surface  of  the 
soil  layer. 
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2.3  Boundary  Conditions 

In  common  with  the  usual  practice  employed  when  analysing, 
for  example,  rolled  filled  layers  of  earth  dams  or  strip  foundations 
or  embankments,  two-dimensional  flow  of  pore  water  is  considered. 

I  he  excess  pore  pressure,  u,  at  the  previous  boundary  surfaces 
(Fig.  2.1)  vanishes  and  at  the  impervious  boundary  surfaces  (an 
insulated  boundary)  the  normal  derivative  of  excess  pore  pressure 
vanishes.  Due  to  symmetry  only  a  half  of  the  cross  section  need 
be  considered. 


MATHEMATICAL  FORMULATION 

2.4  Total  Stress  Distribution 

Changes  in  stresses  take  place  when  a  load  is  applied  to 

the  surface  of  the  soil  mass.  There  is  a  measurable  total  stress 

increase  at  all  points  within  a  "zone  of  influence".  In  the  case 

of  clay  soils,  it  is  of  interest  to  compute  the  instantaneous  excess 

pore  water  pressure  distribution  due  to  the  applied  load.  With 

time,  this  fluid  stress  will  dissipate,  throwing  increasing  amounts 

of  the  applied  stress  onto  the  soil  skeleton  in  the  form  of  effective 

stress.  Based  on  the  assumptions  employed  in  linear  elasticity  and 

referring  to  the  co-ordinate  system  shown  in  Fig.  2.1  the  induced 

stresses  a  ,  a  and  t  acting  at  a  point  may  be  determined  for  a 
X  z  xz 

strip  load  resting  on  the  surface. 

The  total  induced  stresses  at  a  point,  [S]  within  the 
compressible  soil  layer  may  be  separated  into  dilatational  and  de- 
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FIG.  2.1  GEOMETRY  OF  THE  PROBLEM  AND  THE  BOUNDARY  CONDITIONS  CONSIDERED 


viatoric  components  [D‘]  and  [D" ]  respectively  and  written  as 


CD 

O 

O 

_ l 

a  -  0  0  t 

x  xz 

1 — I 

CO 

1 _ 1 

II 

1 — 1 
o 

1 _ 1 

+ 

1 — 1 
o 

1 _ 1 

II 

o  e  o 

+ 

0  eye  0 

CD 

O 

O 

t.,  0  a  -0 

-4 

XZ  2 

o  +o  +o 

where  6  =  (— — | — -)  2.2 

and  Toct  =  J  / (°x-°y)2  +  (°y-°z)2  +  (oz'ax)2  +  6  4  2-3 

For  plane  strain  problems 


T 


xy 


2.4 


and 


a  =  0.5  cos^  cj)' 


(a  +a  )  (Bishop,  1966) 

X  iL 


2.5 


2.5  Initial  Excess  Pore  Pressure  Distribution 

When  a  load  is  applied  to  the  saturated  soil  mass,  the 
entire  load  is  taken  up  by  the  water  at  the  instant  of  loading. 

The  pore  pressure  generates  a  pressure  equal  to  the  applied  load. 

The  process  of  consolidation  depends  upon  the  rate  at  which  excess 
pore  pressures  are  dissipated.  One  very  important  factor  which 
governs  this  dissipation  process,  is  the  initial  excess  pore  pressure 
distribution  in  the  soil. 


- 
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The  theory  of  elasticity  can  be  used  to  obtain  the  initial 
distribution  of  excess  pore  pressure  in  the  clay  soil.  The  pore 
pressure  distribution  introduces  spatial  gradients  in  the  water 
which  give  rise  to  body  forces  in  the  soil. 

For  an  ideal  saturated  soil  with  an  elastic  skeleton, 
the  initial  pore  pressure,  u. ,  is  given  by 


u 


l 


=  0 


2.6 


Soils  dilate  while  undergoing  volume  changes  as  a  result  of  shearing 
strains  without  variation  in  the  hydrostatic  stress.  Henkel  (1960) 
considered  such  an  effect  by  modifying  Equation  2.6  and  it  is 
generally  accepted  that  the  following  equation  is  a  good  approxi- 
mati on. 


0  +  a  t 


oct 


2.7 


where  "a"  is  the  Henkel  pore  pressure  parameter  which  is  related 
to  the  Skempton  (1954)  "A"  parameter  by  the  relation 

/2a  =  A  -  j  2.8 

The  convenience  of  using  0,  which  includes  no  shearing 
effects,  and  t  t,  which  is  independent  of  hydrostatic  stress, 
is  well  demonstrated  in  Equation  2.7. 
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2.6  Governing  Equation 

Adopting  the  empirical  relation  between  void  ratio,  e, 
and  hydrostatic  effective  stress,  0',  (Henkel  and  Sowa,  1963),  as 


e  =  eQ  -  N  log 


e' 


10  0'. 


2.9 


the  coefficient  of  compressibility,  m  ' ,  is  given  by 


m. 


1  3e  _  0.434  N 
(1+e)  90 1  (1+e)  0 


2.10a 


where 


e'i  =  (o'xi  +  °'yi  +  °'zi)/3  • 


2.10b 


initial  effective  stress 


0' 


(o'x  +  c'y  +  a 


'z)/3 


2.10c 


and  N  denotes  the  slope  of  e-log© 1  curve. 

During  the  consolidation  process,  (1+e)  varies  with  time 
far  less  than  0'  hence,  (1+e)  may  be  considered  constant  for  any 
load  increment.  With  this  assumption  Equation  2.10a  becomes 


_B 

0' 


B  = 


0.434  N 


2.1  Od 


where 


a  constant 


2 .  lOe 


- 
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For  normally  consolidated  clay,  Davis  and  Raymond  pointed  out  that 
the  coefficient  of  consolidation,  Cy,  varies  much  less  than  the 
coefficient  of  compressibility,  mv<  They  assumed  Cy  to  be  a 
constant.  It  can  be  therefore  assumed  that 


Y  m 

'w  v 


=  constant 


2.11 


where  k  denotes  the  coefficient  of  permeability 
Yw  denotes  the  unit  weight  of  water 
and  Cy1  denotes  the  coefficient  of  consolidation. 
Using  Darcy's  law 


z  3jj 

\  dz 


2.12 


where  vz  denotes  the  velocity  of  flow  in  the  z-di  recti  on 

kz  denotes  the  coefficient  of  permeability  in  the  z-di recti  on 
i  denotes  the  hydraulic  gradient  in  the  z-direction 
and  u  denotes  excess  pore  pressure 

and  substituting  Equation  2 . 1 Od  in  Equation  2.11  the  resulting 
equation  is 


2.13 


Substituting  Equation  2.13  in  Equation  2.12  and  dif¬ 
ferentiating  the  resulting  equation  with  respect  to  z  yields 


ti 
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3v  ~  B  C  '  . 

_ z  _  3  r  V  3u-| 

3z  3z  L— 0  3zJ 


2.14a 


n  r  i  3  r  1  3U-1 

LV  9z  L0 1  3zJ 


2.14b 


d  p  i  r  1  3  U 

B  CV  ‘-01 


1  30'  3u 


3Z2  ( 0 1 ) 2  92  9Z 


] 


2.14c 


in  the  same  manner  similar  expressions  such  as  Equation  2.14c  can 
be  produced  for  the  x  and  y  directions.  These  are, 


8V 

-  B  Cv' 

32u 

1 

30 

3x 

3x2 

(e1)2 

3x 

a? 

t<< 

tl 

‘  B  C/ 

32u 

1 

30 

3y 

3y2 

(e1)2 

2.15 


dU  j 


2.16 


The  rate  of  water  lost  from  an  element  dxdydz  is  given  by, 


B  C„'  {gV  ( 


3v 

3v 

3v 

<, x  + 
v3x 

3y 

3z 

.2 

.2 

3  u  , 

9  u  + 

3  u 

2 

2 

2 

3x 

3y 

3z 

)  .  1  j-36'  3u  +  36'  3u 


( e 1 ) 


2  l3x  3x  3y  3y 


+  |ii|i]}  dxdydz 


2.17 


The  strain  developed  within  the  soil  element  is  given  by 


(e0-e) 


2.18a 
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where  e  denotes  strain 

and  eQ  denotes  the  void  ratio  corresponding  to  zero  strain  and 
to  the  stress  6. ' . 

Making  use  of  Equation  2.9,  strain  can  be  rewritten  as 

G  =  ( l+eQ)  log10  2,18 

The  rate  of  volume  change  can  be  obtained  by  differentiating  the 
strain  with  respect  to  time 


3c  _  N 

9t  n+e0) 


0.434  39' 

e'  at 


dxdydz 


2.19a 


Considering  the  change  in  void  ratio,  e,  is  too  small 
to  affect  the  total  volume,  (1+e),  it  may  be  assumed  to  be  constant. 
Hence  the  Equation  2.19a  may  be  rewritten,  making  use  of  the 
Equation  2.10e,  as 


dxdydz  =  V  If  dxdydz  2.19c 

Equations  of  continuity  demand  that  the  amount  of  water 
lost  must  be  equal  to  the  reduction  in  volume.  Therefore 

9v  3v  3  v,  3- 

(gT  +  a/  +  aC  >  dxdydz  =  at  dxdydz  2'20 


Equating  Equations  2.17  and  2.19c  to  satisfy  continuity, 


' 
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the  resulting  equation  becomes 


1  30' 


e  1  Bt 


_  c  1  { 


,  rZ  -2  rZ 
1  fd  u  ,  3  u  ,  9  uN 

F  +  7T  +  7Tj 


3x 


3y 


Bz 


(0') 


r 36 1  Bu  ,  36 1  3u  ,  30 1  Bu-n 
l3x  3x  By-  By  b7~  dzj} 


2.21 


Equation  2.21  is  the  general  three-dimensional  consoli¬ 
dation  equation  for  the  assumptions  employed  in  Davis  and  Raymond 
theory. 

The  original  Davis  and  Raymond  one-dimensional  equation 


(no  flow  in  the  x-direction  and  the  y-direction)  may  be  derived 

8u_ 

3y 


2  2 

by  considering  37  =  =  0  and  ~ jr  ~  =  0-  Equation  2.21  then 


Bx 


becomes 


6y 


Bo' 


3 1 


z  =  -  r  {  — 

LV  1  a' 


32u 


1 


Ba' 


z  Bu 


z  Bz 


(a'z) 


Bt  Bz 


]} 


which  is  the  same  as  given  by  Davis  and  Raymond  (1965). 
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For  two-dimensional  consolidation  (flow  in  the  y-direction 

Bu 
By 


2 

is  zero),  it  can  be  considered  that  —  =  0  and  — j  -  0  and  the 

3v  8y 


governing  equation  becomes 


39 1 

Bt 


„  ,  r,B2u  ,  32un  1  ( 30 1  Bu  A  86 1  3un  ^ 
"  LV  U7I  „_2'  "  0'  k3x  Bx  Bz  Bz' s 


Bx 


Bz 


2.23 
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For  a  constant  load  with  time,  (a  +  a  +  a 

x  y 

be  considered  as  constant,  hence 


By  definition 


0  =  6'  +  u 

and  differentiating  Equation  2.25  and  substituting  in 
2.24  yields 


99* 

at 


9u_ 

at 


The  stress  components  a  ,  a  and  t  for  the 

x  z  xz 

medium  (Fig.  2.1),  are  (Harr,  1966) 


°x  -  a  [tan-1  ^  -  tan-1  *§■  .  ] 

x  z  z  (x  +z2-B2)2+4B2z2 

q  rj.  -1  x+B  -1  x-B  ,  2zB(x2-z2-B2) 

a,  =  -  [tan  —  -  tan  —  +  o  \  ■?  9 — L 


Z  TT 


(x2+z2-B2)2  +  4B2z2 


and 


=  4q  r  Bxz _ _ 

Txz  77  (x2+z2-B2)2  +  4B2z2 


z  >  0,  I x I  < 


oc 


,)/3  =  e  can 

2.24 

2.25 

Equati on 

2.26 

soi  1 

2.27 

]  2.28 

2.29 
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where  q  denotes  intensity  of  load 

and  B  denotes  half  of  the  loaded  width. 

Using  the  Equations  2.5,  2.25,  2.26,  2.27,  2.28  and  sub¬ 
stituting  in  Equation  2.23  yields 


9u.  =  r  . 

3t  V 


2  2 
,3  u  ,  3  u\ 

\  9  O  /  “ 


3x 


3z 


{[C^ 


z2+(x+B)2 


z2+(x-B)2 


9 u \  3u -i 
3x'  3xJ 


+  c  r  (x-B) 
'z2*(»B)2 


(x+B) 

z2+(x+B)2 


9u -i  3u-| 
3zJ  3z^ 


2.30 


where  C-j  =  q [1  +  0.5  cos2^1]. 

Equation  2.30  is  the  general  governing  equation  for  con¬ 
solidation  of  a  normally  consolidated  soil  with  two  dimensional 
drainage.  This  equation  is  highly  non-linear  and  could  not  be 
solved  with  the  numerical  methods  adopted  in  this  study. 

Assuming  a  constant  load  with  respect  to  time  as  well  as 
with  space  (as  in  the  case  of  thin  clay  deposits)  and  differentiating 
Equation  2.25  with  respect  to  time  and  space  and  substituting  in 
Equation  2.23,  yields 


3u  _ 
3t 


f  1 


u4  + 

3x 


3  u  ^ 

2' 

3z 


<T  [(f)2  +  (f  )2]} 


2.31 


Equation  2.31  represents  two-dimensional  consolidation  of 
normally  consolidated  soils  under  a  constant  load.  A  solution  to 
Equation  2.31  will  be  unique  depending  on  specified  boundary  and 
initial  conditions. 


’ 


Using  the  substitution 
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a.  e'f-u 

w  -  lo9io  (T*  = 


'10  0 


2.32a 


where 


e'-r  =  ( 


a>xf  +  a  *  yf  +  °'zf 


=  final  effective  stress 


and  differentiating  with  respect  to  x 


2.32b 


3w  _  0.434  3u 

3x  0'  3x 


2.33 


2  2 
3  W  p.  r  o  n  f  1  3  U  ,  1  /3U\2*1 

—  =  -  °-434  LqT-  ^“2  +  - — —  (— )  ] 

3x 


3x2  (0')2  3x‘ 


2.34 


Si mi larly , 


3w  _  0.434  3u 

3z  “  0'  3z 


2.35 


=  .  0.434  [A  ^  +  — Ly  (|^)2] 
3z2  6  3z2  (d') 


2.36 


yields 


Differentiating  Equation  2.32a  with  respect  to  time. 


3w  _  0.434  30 1 
3t  0'  3t 


2.37 
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Substituting  Equations  2.34,  2.36  and  2.37  into  Equation 
2.31  gives  the  simpler  version  of  the  differential  equation  in 
terms  of  the  function,  w 


9w  =  r  , 
at  LV 


[8^w  +  a^w-j 


ax‘ 


az 


2.38 


Equation  2.38  is  identical  to  that  of  the  Terzaghi 
Equation  1.3,  and  can  be  solved  in  the  same  way. 

For  the  problem  considered  herein  (a  strip  load  of  width, 
2B,  of  uniform  load  intensity,  q,  resting  on  a  half-space)  the 
boundary  conditions  are  similar  in  terms  of  u  and  w  (Fig.  2.1). 
These  conditions  are  given  in  Table  1. 
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CHAPTER  III 

ANALYSIS  BY  NUMERICAL  METHODS 


3.1  General 

Closed-form  solutions  of  the  non-linear  equation  for 
specified  boundary  and  initial  conditions  are  difficult  to  obtain 
and  alternate  methods  must  often  be  used,  for  example,  numerical 
methods . 

The  expansion  of  the  function  by  Taylor's  series  provides 
the  fundamental  concept  in  the  finite  difference  approximation  of  a 
differential  equation.  By  representing  the  derivatives  as  finite 
differences,  a  system  of  differential  equations  is  converted  to  a 
system  of  difference  equations  and  the  problem  reduces  to  the  solu¬ 
tion  of  a  set  of  simultaneous  equations. 

Most  numerical  schemes  call  for  step  by  step  procedures. 

In  the  process  of  solution  of  differential  equations,  errors  are 
introduced  at  each  step  of  the  calculation  due  to  the  approximations 
in  finite  difference  and  numerical  calculations.  The  finite  dif¬ 
ference  scheme  must  be  so  formulated  that  the  error  growth  with 
each  step  is  tolerable.  If  error  grows  too  rapidly,  the  scheme  is 


unstable. 


' 
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3.2  Alternating-Direction  Implicit  Method  (ADI) 

The  alternating-direction  implicit  (ADI)  method  was  intro¬ 
duced  by  Peaceman  and  Rachford  (1955)  for  use  in  problems  involving 
a  large  number  of  internal  nodal  points.  In  this  method  the  size 
of  the  matrix  to  be  solved  at  each  time  step  is  reduced  at  the  ex¬ 
pense  of  solving  the  reduced  matrix  many  times  for  each  time  step. 

An  alternating-di recti  on  implicit  procedure  requires  a 
line-by-line  solution  of  small  sets  of  simultaneous  equations.  These 
equations  may  be  solved  by  a  non-iterative  direct  method.  Direct 
methods  for  solving  sets  of  equations  on  a  computer  are  usually 
based  on  Gaussian  elimination. 

The  finite  difference  representation  for  the  two-dimensional 
consolidation  equation, 


2  2 

3u  _  3  u  ^  3  u  -31 

^%z2  3x2 

in  the  ADI  method  (Fig.  3.1)  is  given  by 

{u(L ,M,N+1 )  -  u(L,M,N)}/AT 
=  {u(L-l ,M,N+1 )  -  2u(L,M,N+l)  +  u(L+l ,M,N+1 ) }/(Az)2 
+  {u(L,M-l ,N)  -  2u(L,N,N)  +  u(L,N+l ,N)}/(Ax)2  3.2 


This  is  used  to  advance  the  solution  from  the  N  to  the  (N+l)  time 
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step  and  the  equation. 


{u(L  ,N  ,N+2)  -  u(L  ,N,N+1 ) }/AT 
=  {u (L- 1  ,N,N+1)  -  2u(L,N,N+l)  +  u(L+l ,M,N+1 ) }/(Az)2 
+  {u(L,M-l ,N+2)  -  2u(L,M,N+2)  +  u(L,M+l  ,N+2)}/(Ax)2  3.3 

is  used  to  advance  from  the  (N+l)  to  the(N+2)  time  step. 

This  method  consists  of  replacing  only  one  of  the  second- 
32u 

order  derivatives,  — j  ,  say,  by  an  implicit  difference  approximation 

3z 

in  terms  of  unknown  values  of  u  from  the  (N+l)  time-level,  the 

2 

other  second-order  derivative,  — j  ,  being  replaced  by  an  explicit 

3x 

finite  difference  approximation  in  terms  of  known  values  of  u  from 
the  Nth  time-level. 

The  solution  is  advanced  to  the  (N+2)  time- level  by  re- 
32u 

placing  by  an  implicit  finite-difference  approximation  and 
3x 

32u 

by  an  explicit  one.  The  use  of  each  of  the  above  equations 
3z^ 

at  each  time  step  leads  to  P  sets  of  P  simultaneous  equations. 

Because  these  are  solved  alternately  by  sets  of  rows  and  sets  of 
columns,  the  method  may  be  considered  as  a  line  method  with  alter¬ 
nating  directions. 


1 


3.3  Finite  Difference  Representation 

Using  the  non-dimensional  variables, 


U  = 


_u 

el 

i 


3.4 


Z  = 


2 

H 


3.5 


X  = 


H 


3.6 


and 


T  = 


V  t 

H2 


=  time  factor 


3.7 


Equation  2.31  can  be  rewritten  as 


3U  _  f,32U  32IK  .  1  r,3U>2  I  r3Uv2-,, 

3T  ‘  {(^2  +  +  -  [(3T>  +  {3X>  ]1 


(or1-  -  U) 


9  i 


e  f 


3.8a 


Defining  “r1- =  DR  (stress  increment  ratio)  Equation  3.8a 
w  1 


becomes 


3U  r/9  U  j  3  U\  x  1  r/9U\2,/9U\2n  oou 

FT  =  {(“2  +  ^  +  TdFDT  3Z"  +  ^~3X  ]}  3-8b 


The  governing  Equation  3.8b  can  be  expressed  in  finite 
differences  (making  use  of  the  ADI  method)  as 


(U(L ,M,N+1 )  -  U(L ,M,N) }/AT 
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=  {U(L-1,M,N+1)  -  2U(L,M,N+1 )  +  U(L+1 ,H,N+1 ) }/(AZ)2 
I  +  UJ(L,M-1 ,N)  -  2U(L,M,N)  +  U(L,M+1 ,N)}/(AX)2 
+  { 1  / [DR ( L , M)  -  U(L,M,N)]}  {[(U(L+1  ,M,N)  -  U(L-1 ,M,N) )/2AZ]2 
+  [(U(L,M+1,N)  -  U(L,M-1  ,N))/2AX]2} 
and  for  the  next  time  step,  as 

U(L ,M,N+2)  -  U(L,M,N+1 )  /AT 
=  {U ( L- 1 ,M,N+1 )  -  2U(L,M,N+1)  +  U(L+1 ,M,N+1 ) }/(AZ)2 
+  {U(L,M-1 ,N+2)  -  2U(L ,M,N+2)  +  U(L,M+1 ,N+2)}/(AX)2 
+  {1/[DR(L,M)  -  U(L,M,N+1 )]}{[(U(L+1 ,M,N+1)  -  U(L-1 ,M,N+1 ) )/2AZ]2 
+  [(U(L ,M+1 ,N+1 )  -  U(L,M-1,N+1))/2AX]2}  3.10 


Since  an  impervious  boundary  exists  at  L  =  1  and  M  =  1 
(Fig.  3.1)  there  is  no  hydraulic  gradient;  therefore 


U(0 ,M,N+1 )  =  U(2 ,M,N+1 ) 


3.11 


U(L  ,0  ,N+1 )  =  U(L,2,N+1) 
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At  L  =  I  and  M  =  J  a  pervious  boundary  exists  and  as  such  no  excess 
pore  pressure  is  possible,  therefore 


U(I ,M,N+1 )  =  0 


U(L,J,N+1)  =  0 


3.12 


The  initial  average  effective  stress  is  calculated  from 
Kq  conditions.  The  initial  excess  pore  pressure  u^  due  to  the  ex¬ 
ternal  load  is  calculated  from  Equation  2.7.  The  stress  increment 
ratio  is  given  by 


(initial  average  effective  stress  +  initial  excess  pore  pressure) 

(initial  average  effective  stress) 

The  governing  equation  is  converted  into  a  set  of  simul¬ 
taneous  algebraic  equations  using  finite  differences.  These  equations 
can  be  solved  by  Gaussian  elimination. 

The  same  procedure  can  be  adopted  for  the  linear  Equation 
2.38  with  appropriate  initial  conditions.  Since  both  the  equations 
are  mathematically  equal  they  should  yield  the  same  solution  for  the 
same  problem.  Unfortunately,  the  non-linear  Equation  2.31  could  be 
solved  only  for  a  low  range  of  stress  increment  ratio  (constant 
ratio  up  to  5)  whereas  the  linear  Equation  2.38  could  be  solved  for 
any  range  of  stress  increment  ratio.  Where  comparison  of  stress 
increment  ratios  is  possible  agreement  is  found. 


‘ 
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Experience  indicates  that  Equation  2.31  is  not  stable 
for  the  numerical  procedure  used  herein.  This  is  also  true  for 
the  general  governing  Equation  2.30.  In  order  to  extend  the  range 
of  solution,  the  linear  equivalent  must  be  solved. 
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CHAPTER  IV 

RESULTS  AND  DISCUSSION 

Equation  2.31  represents  the  general  problem  of  two- 
dimensional  pore  pressure  dissipation  of  a  strip  load  of  width  2B, 
of  uniform  load  intensity  q,  (Fig.  2.1)  resting  on  a  half-space. 

The  effective  angle  of  internal  friction  (<J>')  is  considered  as  25° 
and  the  submerged  weight  of  soil  (y1)  as  55  pcf.  For  analysis,  the 
non-linear  Equation  2.31  is  transformed  to  the  linear  equation  2.38 
by  defining  w  =  1  og-j q  0'/O'f.  Equation  2.38  also  yields  the  progress 
of  consolidation  in  two- dimensions  and  it  may  be  solved  for  a  given 
set  of  boundary  conditions. 

Due  to  the  variation  of  stress  increment  ratio  (DR)  within 
the  defined  clay  region  0  <  x  <  3B  and  0  <  z  <  H,  the  excess  pore 
pressure  will  vary  spatially  at  time  t  =  0.  The  initial  condition 
is  given  by  (from  Equation  2.32a) 


w(x,z,0)  =  loglo(0'i(x,z)/0,f(x,z)) 
=  log10(l/DR(x,z)) 


4.1 


The  boundary  conditions  for  the  clay  layer  (only  half 
the  section  need  be  considered  for  the  analysis  because  of 


. 
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symmetry)  at  any  time  t  >  0  are 

3w  n 

3^  =  0  at  x  =  0 

3w  n  , 

97=0  at  z  =  H  4.2 

w  =  0  at  z  =  0  and  x  =  3B 

For  the  given  initial  and  boundary  conditions  (Equations 
4.1  and  4.2  respectively)  Equation  2.38  can  be  solved  as  an  ordinary 
diffusion  equation  using  numerical  finite  difference  methods  such 
as  ADI.  The  average  degree  of  consolidation  U,  along  the  vertical 
section  passing  through  the  centreline  is  given  by 

/  ut dz 

U  =  1 - 5 -  4.3 

/  ^  dz 

where  u^.  denotes  excess  pore  pressure  at  time  t 

and  u.  denotes  initial  excess  pore  pressure. 

\ 

The  average  degree  of  consolidation,  U,  is  calculated 
at  each  time  step  by  integration  making  use  of  Simpson's  rule  for 
the  even  numbers  of  ordinates  in  the  range  of  integration. 

The  boundary  specified  by  x  =  BB  and  0  £  z  £  H  is  studied 
to  determine  its  effect  on  the  computed  results.  The  results  of 
the  relationship  between  average  degree  of  consolidation  U  and  time 
factor  shown  in  Fig.  4.1,  Fig.  4.2  and  Fig.  4.3  indicate  that  a 


. 


Time  Factor 
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AVERAGE  DEGREE  OF  CONSOLIDATION  versus  LOG  TIME  FACTOR 


Time  Factor 
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AVERAGE  DEGREE  OF  CONSOLIDATION  versus  LOG  TIME  FACTOR 
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VALUES  OF  AVE.  DEGREE  OF  CONSOLIDATION  FOR  A=1.0,  q=1000  psf  AND  DIFFERENT  VALUES  OF  3 
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VALUES  OF  AVE.  DEGREE  OF  CONSOLIDATION  FOR  A=1.0,  q=2000  psf  AND  DIFFERENT  VALUES  OF  3 
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distance  4B(  i.e.  3=4)  is  adequate  in  order  to  neglect  the  in¬ 
fluence  of  the  boundary.  It  is  of  interest  to  note  that  Krizek  et 
al  (1969)  have  drawn  similar  conclusions  in  their  study.  Com¬ 
puted  results  are  given  in  Table  2,  Table  3  and  Table  4. 

In  their  one-dimensional  study,  Davis  and  Raymond  (1965) 
noted  that  the  average  degree  of  consolidation,  U,  decreases  with 
the  increase  of  stress  increment  ratio  for  a  given  time  factor,  T  . 
In  the  present  study  the  shape  of  the  curves  in  Fig.  4.4  indicate 
the  same  tendency.  As  the  applied  load  increases,  the  stress  incre¬ 
ment  ratio  also  increases  with  the  result  that  the  average  degree 
of  consolidation  decreases  with  increasing  load.  Computed  results 
in  Table  5  show  that  the  maximum  variation  in  the  average  degree 
of  consolidation  for  loads  of  500  psf  and  1000  psf  is  about  4.5% 
and  the  same  variation  is  found  for  the  loads  1000  -  2000  psf  and 
2000  -  4000  psf. 

In  this  study  it  is  assumed  that  Equation  2.7,  as  sug¬ 
gested  by  Henkel  (1960),  is  a  good  approximation  for  representing 
the  initial  excess  pore  pressure  u^ ; 


For  an  ideal  saturated  soil  with  an  elastic  skeleton,  the 
initial  excess  pore  pressure  is  given  by  Equation  2.6 

u.j  =  0 


.fiq  Q0O>  -  OOOS 


Time  Factor 
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FIG  4.4  AVERAGE  DEGREE  OF  CONSOLIDATION  versus  LOG  TIME  FACTOR 


TABLE  5 
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To  determine  the  effect  on  the  relationship  between 
average  degree  of  consolidation  u  and  time  factor  by  using  Equation 
2.7  rather  than  Equation  2.6,  calculations  have  been  performed  for 
the  problem  shown  in  Fig.  2.1. 

It  is  found  that  the  Henkel  pore  pressure  parameter  "a" 
has  an  insignificant  effect  on  the  rate  of  consolidation.  This  is 
clearly  demonstrated  by  Fig.  4.5,  Fig.  4.6  and  Fig.  4.7.  It  appears 
reasonable  to  conclude  that  modifications  to  the  general  approach 
to  take  into  account  non-elastic  effects  on  the  initial  pore  pressures 
in  terms  of  rate  of  consolidation  are  not  worthwhile  for  practical 
purposes.  The  possible  reason  for  this  may  be  that  the  non-elastic 
effects  constitute  only  about  25%  of  the  elastic  effects  on  the 
distribution  of  initial  excess  pore  pressure.  It  has  already  been 
noted  that  the  maximum  variation  in  the  average  degree  of  consoli¬ 
dation  for  the  range  of  loads  employed  is  about  4.5%.  The  variation 
in  average  degree  of  consolidation  is  about  2%  for  an  extreme 
range  of  "a"  (Table  6) . 

Although  the  influence  of  non-elastic  effects  on  the  rate 
of  consolidation  in  insignificant,  however,  the  amount  of  settlement 
will  be  influenced. 

It  is  known  that  the  intensity  of  stress  in  the  soil  medium 
due  to  external  load,  decreases  with  increase  of  depth  and  with  de¬ 
crease  of  load.  It  controls  the  variation  of  initial  excess  pore 
pressure  within  the  medium.  Since  the  initial  excess  pore  pressure 
at  the  mid-point  of  the  bottom  boundary  is  less,  flow  takes  place 
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toward  the  mid-point  due  to  the  difference  in  hydraulic  gradient. 
This  process  increases  the  pore  pressure  from  the  initial  value  and 
decays  as  time  progresses.  The  effect  is  illustrated  in  Fig.  4.8. 
It  may  be  noted  that  the  increase  of  pore  pressure  is  less  as  load 
increases.  For  a  load  of  500  psf  the  increased  pore  pressure  is 
about  4 %  while  for  4000  psf  it  is  less  than  1%. 

Even  though  the  increase  of  pore  pressure  is  within  a 
tolerable  range  from  the  point  of  view  of  safety,  it  is  of  interest 
to  note  that  it  is  possible  for  the  phenomenon  to  occur.  Some 
results  are  synthesized  in  Table  7.  Figure  4.9,  Fig.  4.10,  Fig. 
4.11  and  Fig.  4.12  display  lines  of  equal  pore  pressure  at  dif¬ 
ferent  stages  of  consolidation. 
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MAXIMUM  U/Ui  versus  LOG  TIME  FACTOR 


VALUES  OF  MAXIMUM  U/Ui  FOR  A=1.0,  3=1.0  AND  DIFFERENT  VALUES  OF  q 
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LINES  OF  EQUAL  PORE  PRESSURE  AT  TIME  FACT0R=0.0  (U%=0.0) 


q= 1000  psf 
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FIG.  4.10  LIMES  OF  EQUAL  PORE  PRESSURE  AT  TIME  FACTORS. 2  (U%=30.523) 
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FIG.  4.11  LINES  OF  EQUAL  PORE  PRESSURE  AT  TIME  FACTORS. 4  (U%=61 .336) 


q=1000  psf 
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FIG.  4.12  LINES  OF  EQUAL  PORE  PRESSURE  AT  TIME  FACTORS. 9  (U%=95.719) 
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CHAPTER  V 
CONCLUSIONS 

On  the  basis  of  a  limited  computer  study  on  the  two-dimensional 
dissipation  of  pore  pressures  in  a  normally  consolidated  soil  under 
a  strip  load  of  uniform  load  intensity,  the  following  conclusions 
may  be  drawn: 

1)  The  rate  of  dissipation  of  excess  pore  pressure  is  al¬ 
most  independent  of  pore  pressures  induced  by  shear 
stresses. 

2)  Average  degree  of  consolidation,  F,  decreases  with  in¬ 
crease  of  load  intensity. 

3)  Increase  of  initial  excess  pore  pressure  is  possible 
at  the  mid-plane  after  consolidation  commences. 

The  conclusions  presented  herein,  are,  for  the  assumption 
that  the  total  stress  remains  constant  with  respect  to  time  as  well 
as  with  space  which  is  not  true  for  the  general  case  (constant  only 

with  time).  As  Schiffman  et  al  (1969)  point  out  " . numerical 

techniques  for  solving  non-linear  problems  are  not  well  defined", 
regrettably  the  general  case  could  not  be  solved  with  the  numerical 
techniques  adopted  in  this  study. 
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COMPUTER  PROGRAM  AND  USAGE 
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APPENDIX 

The  computer  program  presented  herein  describes  the 
process  ot  consolidation  in  two-dimensions  (for  thin  layers)  as 
given  by  the  Davis  and  Raymond  theory.  The  program  calculates  the 
ratio  of  pore  pressure  to  the  initial  pore  pressure  at  all  nodal 
points  and  the  average  degree  of  consolidation  at  each  time  step. 
The  program  is  in  FORTRAN  IV  language  and  may  be  used  directly  on 
computers  of  Type  360/67. 

The  equation  governing  the  dissipation  of  excess  pore 
pressure  in  two-dimensions  is  given  by  Equation  2.38.  The  method 
of  analysis  by  the  ADI  technique  is  described  in  Chapter  III. 

The  program  is  applicable  for  a  strip  load  of  width  2B , 
of  uniform  load  intensity  q,  imposed  on  the  surface  of  a  half-space. 
The  program  takes  into  consideration  the  stress  increment  ratio 
(DR)  at  all  nodal  points.  The  distribution  of  initial  excess  pore 
pressure  depends  on  the  stress  increment  ratio. 

The  entire  program  consists  of  a  main  and  two  subroutines. 
The  main  program  reads  the  input  data  and  sets  up  the  initial  and 
boundary  conditions.  The  main  program  also  regulates  the  channel 
to  the  subroutines  by  setting  up  appropriate  conditional  statements. 

Subroutine  JAM  computes  the  stress  components  due  to  the 
external  loading  and  calculates  the  stress  increment  ratio.  Sub¬ 
routine  LEENA  gives  the  ratio  of  pore  pressure  to  initial  pore 
pressure  at  all  nodal  points. 

For  the  substitution  employed  (Equation  2.32a)  the  initial 
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condition  at  time  zero  is  (i.e.  N  =  1) 


A-2 


e'_. 


w(L,M,l)  log]0  qi  -  1  °g-j q  [1/DR(L ,M) ] 


A.  1 


For  any  time,  say  N  (Fig.  3.1) 


w(L,M,N)  =  login  { 


0 1  f  -  u(L,M,N) 


10  L  6 


A. 2 


Combining  Equation  A. 2  and  Equation  2.25,  yields 


u(L,M,N) 

Uj (L ,M) 


=  DR(L,M)  (1 


e2.303  w(L ,M,N) 


}/[DR(L,M)  -  1]  A. 3 


Equation  A. 3  is  the  ratio  of  pore  pressure  to  the  initial 
excess  pore  pressure.  Pore  pressure  at  any  node  can  be  calculated 
by  multiplying  Equation  A. 3  with  the  corresponding  initial  excess 
pore  pressure. 

The  average  degree  of  consolidation  on  the  basis  of  the 
centre-line  integral  method  is  calculated  making  use  of  Simpson's 
rule.  The  calculated  values  are  printed  out  according  to  the  out 
put  formats. 


Computer  Program  Usage 
Input  Data 

To  calculate  the  stress  increment  ratio  and  assign  the 
boundary  and  initial  conditions,  the  input  data  cards  read  the  in¬ 
tensity  of  load,  the  soil  parameters,  space  steps  in  both  the  di- 


recti ons  and  the  time  factor  steps. 

The  READ  variables  are: 

GAMA,  Q,  DT,  AB,  I,  J,  K  and  DZ ,  DX. 

The  FORMAT  is, 

4F10.5,  3110  and  2F10.5. 

Output  Information 

The  following  information  is  computed  and  printed  by  th 
program  at  required  time  intervals: 

1)  the  ratio  of  pore  pressure  to  the  initial  excess  pore 
pressure  at  all  nodal  points, 

2)  the  time  factor  and  the  average  degree  of  consolidation 
(one-dimensional  integral)  along  the  centre-line  and 

3)  values  of  pore  pressure  at  the  required  node. 


' 


. 

'  - 


o  o  o  on 


C 

C 

C  DAVIS  AND  RAYMOND  TWO-DIMENSIONAL  ANALYSIS 

c 

c 

C  THIS  PROGRAM  CALCULATES  THE  RATIO  OF  EXCESS  POPE  PRESSURE  TO  THE 

C  INITIAL  EXCESS  PORE  PRESSURE,  PORE  PRESSURE  AND  THE  AVEREGE  DEGREE 

C  OF  CONSOLIDATION  BY  ONE- DI ME NS  I ONAL  INTEGRAL  METHOD 

C 

c 

c 

C  ANALYSIS  FOR  A  STRIP  LOAD  OF  WIDTH  2B  OF  UNIFORM  INTENSITY, Q  RESTING 

C  ON  THE  HALF-PLANE 

C 

C  BA=HALF  THE  WIDTH  OF  LOADING 

C  Q=  FOOTING  LOAD  IN  LBS. /S. FT 

C  SUBMERGED  WEIGHT  OF  SOIL  (GAMA)  =55.0  PCF 

C  PI  =  EFFECT  IV  E  ANGLE  OF  INTERNAL  FRICTION 

C  H=DEPTH  OF  THE  SECTION 

C  AB=  SKEM  PTON  PORE  PR.  PARAMETER - A - 

£  $***$*****#***#*=)<*#*#***<«*  ***$#£  **#**$$  «**««*«*«  $*$*#*** 
c 
c 
c 

C  TREATMENT  IN  A  DIMENSIONLESS  FORM 

C 

C  THIS  IS  DOME  WITH  RESPECT  TO  DEPTH, H  SO  THAT  THE  VARIATION 

C  IN  THE  QUANTITY  IN  THE  VERTICAL  DIRECTION  RANGES  FROM  0.0  TO  1.0 

C  AND  IN  THE  HORIZONTAL  DIRECTION  RANGES  FROM  0.0  TO  BETA 

C 

C  TIME  FACTOR  IS  DIMENSIONLESS  AND  IS  EQUAL  TO  CV* *T I  ME /H**2 

C 

c 

c 

COMMON  U (11,41, 21)»DR(11, 41), A(11),D(11)»T(21),B(41),E(41), 

CGAMA , Q , BA, AB , DZ , DX , J , IA,JA,AZ(41),UZ 
COMMON  /  SQ  /  I U  J (  11,41) 


1000  RE AD (5, 101, END=999 )  GAM A , Q, DT , AB , I , J , K 
101  FORMAT (4 F1C. 5 ,3  1 10  ) 

READ(5,201)  DZ, DX 
201  FORMAT ( 2  F 1 0 • 5 ) 

WR ITE(6,541 )  GAMA, Q, DT ,AB, I , J ,K 

541  F0PMAT(4X, *  GA  MA  =  '  ,F10.5,2X, 'Q='  ,F10.5,2X, ' DT  =  ' 
CF10.5,2X,'  1='  ,  1 10 , 2X , '  J  = ' ,I10,2X,'K=', 110) 

WRITE (6, 542)  DZ,DX 

542  F0RMAT(4X,,DZ=,,F10.5,2X,*DX=,,F10.5) 


F10.5,2X,  ' AB=  * , 


INITIAL  AND  BOUNDRY  CONDITIONS 


B  A=  1  •  0 


- 


1  ' 


non  ooo  ooo 


BET  A= ( J- 1 )*DX/BA 

I  A  =  I  —  1 

J  A= J- 1 

DO  10  N= 1 »  K 

DO  10  L= 1 » I 

DO  10  M= 1 ,  J 

IF(N.EO.l)  GO  TO  30 

IFIM.EQ. J.OR.L.EQ.I)  U(L»M,N)=0.0 

GO  TO  10 

30  I F ( L • EQ.  I  )  GO  TO  31 
CALL  JAM(L,M) 

U(L,  V.f  N)  =0.434*  A  LOG  (  (  1 .  ) /DR < L , M )  ) 

GO  TO  10 

31  U(L,M,N)=0.0 
10  CONTINUE 

PHYSICAL  CONSTANTS 

RA=  DT / ( DZ )**2 
R  3= DT / ( DX)**2 
WRI T  E ( 6  » 1 10 )  RA,RBtBETA 

11C  FORMAT ( // /4X  ,  » RA  =  *  , F 1 0. 4, 2X ,  * RB= • , F 10 . 4, 2X , • B ET A= • , F 10 . 4 ) 
TWO  DIMENSIONAL  DISSIPATION  BY  A  D  I  MFTHOD 
KK=K- 2 

DO  121  MM= 1 t 15 
DO  11  N= 1 » KK »  2 

DISSIPATION  IN  THE  VERTICAL  DIRECTION 

Ad)=2.*RA/( 1  •  +2 • *R  A  ) 

A A=  (  1 .  -2  .  *R B  )  *11  (  1, 1  ,N)+2.*RB*U(  1,2*N) 

BB=1.+2.*RA 
D ( 1 ) =AA/BB 
DO  12  L  =  2  » I  A 
M=1 

A  (  L  )  =  R  A/  { 1 . +2.*RA-RA*A( L-l ) ) 

A  A= ( l.-2.*RB)*U(L,M,N)+2.*RB*U(L ,M+1?N) 

BB=1 .+2.*RA-RA*A (L-l ) 

D(L )=(RA*D(L-l)+AA) /BB 

12  CONTINUE 

U( IA, ltN+1) =D( I  A) 

DO  13  L=2,IA 

I  S= I  —  L 

M=1 

U(  I  S»M»N-d  )=D(  IS  )  +A(  I  S)  *U<  I  S  +  l  »M  ,N+  1) 

13  CONTINUE 

DO  14  M=2 »  JA 
L  =  1 

A(1  )  =2.*PA/ (1.^2  .  *P  A ) 

AA=RB*U (L  »M  -  1  »  N  )  *-  ( l.-2.*RB)*U(LtM»N)  +  RB*U(L»M+l »N) 
BB=1.+2.*PA 
0(1) =AA/BB 
DO  15  L  =2  » I  A 


* 


■ 


O  O  o  o 


A(L ) =RA/ ( 1 .  +  2 • *R A-RA$A ( L-l ) ) 

AA=RB*U ( L  »  M- 1 , N) +  d .-2 .*RB  >  *U(L , M,N )+RB*U( L ,  M+l ,  N  J 
BB=1.+2.*PA-RA*A(L-1  ) 

D(L)  =  (P.A*D(L-1  )+AA)/BB 
15  CONTINUE 

U( I A ,M,N+1 ) =D( I  A) 

DO  16  L=  2  »  I  A 
I S= I - L 

U( I S»  MfN+1 )  =  D( IS ) +A ( I S) *U( IS  +  lfMfN+1) 

16  CONTINUE 
14  CONTINUE 

C 

C  DISSIPATION  IN  THE  HORIZONTAL  DIRECTION 

C 

B  (  l)=2.*RB/(  l.+-2.*R8) 

AA= ( 1.-2 ,*R A ) *U( It 1»  N  +  1 )+2.*RA*U(2tl»N+l) 

BB=1.*2.*RB 
E(1 ) = AA/BB 
DO  17  M=2  »  J A 
L  =1 

B(M)=RB/d.  +  2.*RB-RB*B(M-l)  ) 

A A=  (  1.-2  .*R A ) *U ( L  »  M , N  + 1  I +2 . *RA*U ( L+ 1 , M , N  +  l ) 

BB= 1 . +2.*RB-RB*B ( M— 1 ) 

F (M )=( R8*F ( M- 1 ) +AA) /BB 

17  CONTINUE 

U(l, JA,N+2)=E( JA) 

DO  18  M  =  2  »  J A 
L  =  1 

I  P  =  J  —  M 

U(L, IP,N+2)=E( IP)+B( IP)*U(L, IP+l,N+2) 

18  CONTINUE 

DO  19  L  =  2  t I  A 
M  =  1 

B  (  1  )=2.*RB/  (1.^-2  .*RB) 

AA=RA*U(  L-l, M.N+l ) +(1 .-2 . *RA) *U( L, M,N+1  )  *R  A*U  (  L+  1 ,  M  ,  N  + 1 ) 

BB=1 . +2.*RB 
E(l)= AA/BB 
DO  2  2  M=  2 , J  A 

B(M)=RR/( l.+2.*RB-R8*B(M-l) ) 

AA  =  R  A*  U ( L-l , M, N  +  l )  +  { 1 .-2. *RA) *U( L, M,N+1  )  +RA*U { L  +  1 , M , N  + 1 ) 

BB= 1 .+2.*RB-RB*B ( M— 1 ) 

E(M)=(RB*F(M-1 ) +  AA ) / B  3 
22  CONTINUE 

U(L , JA ,N+2 )=E ( JA ) 

DO  2  3  M= 2  »  J  A 
I P=  J-M 

U(L,  IP,N*2)  =  E(  IP)+B(  IP)*U(L,  IP-HfN  +  2) 

23  CONTINUE 

19  CONTINUE 
11  CONTINUE 

CALCULATION  OF  THE  RATIO  OF  PORE  PRESSURE  TO  THE  INITIAL  PORE  PRESSURE 
AT  ALL  NODAL  POINTS 


DO  24  N= 1 , K 


■ 


. 


.. 


ooooooo  oooooo 


DO  24  L= 1 » I  A 
DO  24  M= 1 ,  J 
FE=U(L ,M,N)/( 0.434) 

F=  E  X  P  (  E  E  ) 

U I  L  ,  M*  N)  =  (  DR  (  L  »  M  )  *  (  1  •— F  )  )  /  I  DR  (  L  »  M  )  r*  1  •  0  ) 

24  CONTINUE 

DO  2  6  M=  1 ,  J 
IFIMM.GT.l )  GO  TO  27 
UUiM,l)«1.0 
GO  TO  26 

27  U( I *M, 11=0.0 
26  CONTINUE 

$****$*$*  $**********#**,(.£**££$£££ 

OUT  PUT  PRINT 

**********************  *$****:(;#  *$***$*, 


WRITE(6,541)  GAMA,Q,DT,AB,I ,J,K 
DO  8  5  N=  1  *  K , 2  0 
T(N )=DT*( N-l  ) 

WRITE! 6, 102)  TIN) 

102  FORMAT I///50X,' TIME  F ACTOR- » , F9 . 4 ) 

WRI T  E ( 6 , 103  ) 

103  FOPMATI//5X, 'DISTANCE* ,40X, 'RATIO  OF  THE  PORE  PRESSURE  TO  THE 
CINITIAL  PORE  PRESSURE') 

WRITEI6, 105) 

10  5  FORMAT! /7X, » X /Y '  , 8X , ' 0 . 0* , 7X , ' 0 . 1 ■ , 7X , ' 0 . 2 ' , 7X , ' 0. 3  *  ,  7X , ' 0. 4 ' , 7X , 
C  *  C • 5  5  »7X,  'C.o!»7X»  *  0 • 7  *  »7X,  'O.S' »7X, '0.9' »7X»  9l.G'/) 

CALL  LEFNA! 1  ,11 ,  I , N) 

IFIJ.LT .12)  GO  TO  85 
WRITE! 6,306) 

306  FORMAT! /2fiX,'l.l',7X,'1.2',7X,'1.3',7X,'1.4*,7X,'1.5', 

C7X, '1.6'  »7X,  ' 1 . 7 • , 7  X , ' 1.8' , 7X, '1.9'  ,7X,  '2.0'/) 

CALL  LEENA! 12,21 , I, N) 

IF(J.LT.22)  GO  TO  85 
WRITE! 6,3C7) 

307  FORMAT ( /28 X , ’ 2 . 1 • , 7X , • 2 . 2 • , 7X , * 2 . 3 • , 7X , • 2 . 4 ' , 7X , • 2 . 5 • , 

C7X, '2.6' ,7X, ' 2 • 7 ' , 7X, '2.8' ,7X,'2.9' ,7X, '3.0'/) 

CALL  LEENA! 22,31, I, N) 

IF(J.LT .32)  GO  TO  85 
WRI T  E( 6, 308 ) 

3  08  FORMAT! /2  8X, ' 3. 1 ' , 7X , • 3 .2 ' , 7X, • 3.3 ' , 7X , • 3.4* , 7X,  • 3. 5 ' , 

C7X,  '3.6' ,7X,  • 3.7' , 7X, *3. 8* , 7X, '3.9* ,7X, '4.0'  / ) 

CALL  LEENA(32,4l, I , N ) 

85  CONTINUE 

WRITE(6,92) 

92  FORMAT!'  1') 


AVE.  DEGREE  OF  CONS.  FROM  CENTRELINE  INTEGRAL  METHOD 
USING  SIMPSON'S  RULE 

NO  MODIFICATION  IS  REQUIRED  SINCE  THE  NUMBER  OF  ORDINATES 
IN  CONSIDERATION  IS  EVEN 


\ 


' 


O  O  o  o  o  o  o 


WRITE! 6,206) 

206  FORMAT I///50X, • TIME  F ACTOR 1 0 X ,  • AV E . D EGREE  OF  CONSOLIDATION1) 
DO  9  3  N=  1 »  K, 2 
IFIMM.GT.l )  GO  TO  90 
TIN )=DT*( N-l ) 

M=  1 

AZ(N)=U(  1»M»N)-U(I,M»N) 

DO  106  L  =  2  » I  A  »  2 

AZ(N)=AZ(  N)+4.*IJ(L,M,  NI  +  2.<!U(LHt  M,N) 

1C6  CONTINUE 

AZ(N)=AZ(N)*DZ/3. 

AZZ=AZ ( 1 ) 

GO  TO  107 

90  T(N)=DT*(  (N-l  )  +- (  MM- 1  )*(K-1  )  ) 

M=1 

AZIN) =U ( 1,M,N)-U( I »  M »  N ) 

DO  96  L  =  2  » I  A  »  2 

AZ(N)  =  AZ(  N>  +4.*U(L, M, N)+2.*U(L+1  ,Mf  N) 

96  CONTINUE 

A  Z  I N ) =  A  Z I N) *DZ/3. 

107  UZ=1«— AZIN)  /  ML 

UZ=UZ*100.0 
WR I TE I  6 , 93  )  T l N ) » UZ 
9  3  FORMAT  I  / /50X  ,  FI  0.4 ,2CX,F 10.3) 

98  CONTINUE 

CALCULATION  OF  EXCESS  PORE  PRESSURE  AT  THE  MID-PLANE  OF 
BOTTOM  SECTION 

WRITE! 6,799) 

799  FORMAT  I  /  /.3QX  ,  •  T  I  ME  F  ACTOR  *  ,  10X,  •  RAT  I  0  U/UIINITIAL)  'tlOXt'PORE 
CPRESSUP. E  I  PSF  )  *  ) 

DO  80C  N=1,K,20 
IF(MM.GT.l)  GO  TO  900 
Z  T=DT *  I N-l ) 

GO  TO  901 

9C  0  ZT=DT*( (N-l  )  +  I  MM— 1 )*(K-1  )  ) 

901  UP=U( 1 ,1 ,N)*UJ( 1 , 1 ) 

WR ITEI 6, 801 )  ZT, U( 1 , 1 ,N) ,UP 
801  FCRMAT(//30X,F7.3,17X,F9.3,20X,F1C#3) 

800  CONTINUE 


DO  60  L= 1 »  I  A 
DO  60  M  = 1 , J 

AA=( ( DR  I L,M)-1.G )*U(L,M,K)  ) /DR(L,M) 
U(L,M, 1 ) =0. A  34* A LOG! 1 .-A  A) 

60  CONTINUE 
121  CONTINUE 

WRITE  I  6,  305 ) 

805  FORMAT  I  *  1 '  ) 

GO  TO  1000 
999  STOP 
END 


■ 

- 


r>  o  o  o  o  o 


SUBROUTINE  J  AM ( L ,  M ) 

COMMON  U(  11 ,41 ,21) ,DR ( 11,41  ), A( 11 ) , D( 11 ) ,T(21) ,B(41) , E(41 ), 
CGAMA,Q, BA, AB, DZ , DX,J,IA,JA,AZ(41),UZ 
COMMON/ SQ/ I ,UJ< 11,41) 

CONSTANTS  —  ANGLE  (PIM  =  25  DEG  — 

PI=3. 1416*25. /180. 

CO=l.-SIN( PI ) 

ALPA=l.+Q.5*(COS(PI )**2) 

ZETA=0.5*COS( PI )**2 

STRIP  LOADING  - ELASTICITY  APPROACH  — 

UI=( ( l.+2.*CO)*GAMA*( I-L )*DZ)/(3.0) 

AA=AT AN(  (  ( M-l )*DX+BA ) /(  (  I-L ) *DZ )  ) 

BB=  AT  AN  (  ( <M-1)*DX-BA  )/{  (  I-L  )*DZ)  ) 

CC  =  2 • *BA* ( I - L  )  * D Z * (  (  (  M-  l  )  *D X  )  ** 2  -  ( (  I  - L  )  *DZ ) **2-BA*BA) 

DD= (  (  (M-l )*DX)**2+(  ( I-L) *DZ ) **2-R A* B A ) ** 2  +  4. *BA* BA* ( ( (  I-L  )*DZ)**2) 
EE^CC/DD 

S  Z=Q* ( AA-R3-EE ) /3 .1416 
SX=0*( AA-BB+EE)/ 3.141 6 

TX  Z= ( 4. *Q*BA* { M-l ) *  DX  * ( (  (I-L  )  *D  Z ) **  2 ) ) / ( 3 . 14 16*DD ) 

FF=  SORT (  ( SX-ZFTA* ( SZ+SX)  >**2  + (ZETA* ( SX+SZ)-SZ )**2+( SZ-SX ) **2  + 
C6.*TXZ*TXZ ) 

UJ(L,M)=(ALPA*(SZ+SX)/3.  )  +  ( ( AB-(  1  ./3.  )  )*FF/{ SORT ( 2. ) *3. )  ) 

DR ( L  »  M )  =  ( U I +U  J ( L , M )  )/UI 

RETURN 

END 


SUBROUTINE  LEENA ( J1 , J2 , I , N ) 

COMMON  U(  11,41, 21)»DP(11, 41), A(ll),D(ll)»T(21)«B(41)*E(41)f 
C GAM A , Q,BA,AB,DZ,DX,J, IA,JA,AZ(41),UZ 
DO  80  L  =  1 ,  I 
I  L=  I  -  L 
AL  =  D7 *  I L 
I  L=  1  L  +  1 

I F ( J 1  . GT . 1C )  GO  TO  272 

WRITE!  6,  166)  AL  ,  (  UU  L,  M,  N)  ,  M=J  1 ,  J2) 

166  FORMAT (5X, F6.2,5X,11(F6.2,4X)) 

GO  TO  80 

272  WR I TE ( 6, 107 )  AL , ( U ( I L , M, N ) , M=J 1 , J2 ) 

1 C  7  FO  RM  AT  (  5  X  ,  E6 . 2 , 1  5X  ,  1 C  (  F  6 . 2 , 4  X  )  ) 

80  CONTINUE 
RETURN 
END 


